Using OLS Regression to Predict Median House Values in Philadelphia

Introduction

In recent decades, Philadelphia has experienced substantial changes in population dynamics, economic fluctuations, and redevelopment pressures. Among these changes, housing values have become a key indicator of the socio-economic conditions within the city’s neighborhoods. Rising housing values typically reflect a higher presence of wealthier residents or a neighborhood undergoing gentrification, while declining housing values often signal increased economic pressure or decline. Therefore, accurately forecasting the median housing values across Philadelphia is crucial for urban planners and policymakers, helping them to proactively address potential risks such as displacement and disinvestment. Supporting Philadelphia’s long-term vision of fostering vibrant, diverse, and resilient communities is also essential.

To better understand the potential factors influencing housing values in Philadelphia, we identify several key variables: educational attainment, vacancy rates, the proportion of single-family homes, and poverty levels. These factors are closely related to housing market trends and provide insights into the neighborhood’s socio-economic status.

Educational attainment is closely linked to the socio-economic characteristics of neighborhoods. Individuals with higher educational attainment typically earn higher incomes and contribute more to local economies. A higher concentration of individuals with advanced educational backgrounds tends to increase the demand for well-designed housing in affluent neighborhoods. With a consistent supply, housing prices and values are likely to rise. Similarly, a high poverty rate also indicates a declining community. Since many residents cannot afford luxury housing, the housing values in those neighborhoods are likely to be Lower.

High vacancy rates often correlate with declining neighborhoods and reduced median house values. Research shows that vacant properties affect multiple facets of community life, such as housing and neighborhood vitality, crime prevention efforts, and the well-being of commercial districts. As a result, areas with numerous vacant housing units typically see diminished median housing values.

The proportion of single-family homes in an area influences housing values. These homes are generally private and comfortable, making them more desirable in many U.S. housing markets. However, they are relatively common in suburban neighborhoods and may suffer from low accessibility and limited infrastructure, which can negatively affect their property values as well.

In this study, we utilize ordinary least squares (OLS) regression to analyze the relationship between these socioeconomic factors and median house values in Philadelphia. By examining these relationships, we aim to identify critical predictors of median housing values throughout Philadelphia and offer insights for decision-makers and community initiatives.

Methods

Data Cleaning

To predict median house values in Philadelphia, we obtained the original dataset from the United States Census data. The dataset represents census block groups from the year 2000 and initially contained 1,816 observations. The key variables included:

• POLY_ID – Census Block Group ID
• MEDHVAL – Median value of all owner-occupied housing units
• PCBACHMORE – Proportion of residents in the block group with at least a bachelor’s degree
• PCTVACANT – Proportion of housing units that are vacant
• PCTSINGLES – Percentage of housing units that are detached single-family houses
• NBELPOV100 – Number of households with incomes below 100% of the poverty level
• MEDHHINC – Median household income

To refine the dataset for modeling purposes, we applied the following filtering criteria:

1. Retained block groups with a population greater than 40
2. Included only block groups that contained housing units
3. Excluded records where the median house value was below $10,000

Additionally, we removed a specific block group in North Philadelphia that exhibited inconsistencies, with an unusually high median house value (over $800,000) despite a very low median household income (less than $8,000).

After data cleaning, the final dataset contained 1,720 observations.

Exploratory Data Anaylsis

Summary Statistics

We will first examine the summary statistics ( mean and standard deviation (SD)) of key variables in the dataset, including the dependent variables MEDHVAL (Median House Value), and predictors NBELPOV100 (Households Living in Poverty), PCTBACHMOR (% of Individuals with Bachelor’s Degrees or Higher), PCTVACANT(% of Vacant Houses), PCTSINGLES( % of Single House Units).

The mean (\(\bar{X}\)) represents the average value of a variable and is calculated as:

\[ \bar{X} = \frac{1}{n} \sum_{i=1}^{n} X_i \]

where:

• \(X_i\) represents each individual observation
  
• \(n\) is the total number of observations

The mean gives us a single representative value of the dataset. To measure variability, we use the standard deviation (SD), which quantifies how much the values in a dataset deviate from the mean. The formula for the sample standard deviation (\(s\)) is:

\[ s = \sqrt{\frac{1}{n-1} \sum_{i=1}^{n} \left( X_i - \bar{X} \right)^2} \]

where:

• \(X_i\) represents each individual observation
  
• \(\bar{X}\) is the mean of the observations
  
• \(n\) is the total number of observations

A larger standard deviation indicates that the data points are more spread out, while a smaller standard deviation suggests that the data points are closer to the mean.

Distributions

We will also examine the histograms and apply log transformations for key variables to assess whether the transformed variables follow a more normal-like distribution.

Histograms provide a visual representation of how a variable’s values are distributed, helping to identify whether the data follows a normal distribution, is right-skewed, or left-skewed. Linear regression model assume that variables are approximately normally distributed.

• X-axis: The values of the variable (e.g., house prices, income levels).

• Y-axis: The frequency of observations within each bin.

If the histogram is right-skewed, it suggests that a small number of observations have significantly higher values compared to the rest. For variables following a right-skewed distribution, we will apply a log transformation to these variables to improve normality. Since the log transformation is undefined for zero or negative values, we must first check whether any variable contains zero.

• If no zeros presents, we apply the standard log transformation:
  \[ X' = \log_{10} (X) \]
  
• If zeros present, we adjust by adding 1 before taking the log:
  \[ X' = \log_{10} (X + 1) \]
  This ensures that all values remain positive and avoids undefined values.

By comparing the original histograms with the log-transformed histograms, we will assess whether the transformation improves the suitability of the data for predictive modeling.

Correlations

We will analyze correlations between predictors, to detect potential multicollinearity before proceeding with regression analysis, which can distort model interpretations.

Multicollinearity occurs when predictors are highly correlated, which can lead to unstable regression coefficients. It also inflates standard errors, reducing the statistical significance of predictors, and increases the risk of overfitting, as redundant variables do not contribute new information to the model.

The correlation coefficient \(r\) is calculated as:

\[ r = \frac{\sum_{i=1}^{n} (x_i - \bar{x})(y_i - \bar{y})}{\sqrt{\sum_{i=1}^{n} (x_i - \bar{x})^2 \sum_{i=1}^{n} (y_i - \bar{y})^2}} \] In a more concise way, this above formula is also equivalent to the following:

\[ r = \frac{1}{n-1} \sum_{i=1}^{n} \left( \frac{x_i - \bar{x}}{S_x} \right) \left( \frac{y_i - \bar{y}}{S_y} \right) \] where:
- \(X_i\) and \(Y_i\) are individual data points for variables \(X\) and \(Y\), respectively.
- \(\bar{X}\) and \(\bar{Y}\) are the mean values of \(X\) and \(Y\).
- The numerator represents the covariance between \(X\) and \(Y\), while the denominator normalizes the values.

The correlation coefficient \(r\) ranges from -1 to 1: A value of +1 indicates a perfect positive linear relationship, while a value of -1 indicates a perfect negative linear relationship. A value of 0 suggests no linear relationship between the variables, meaning changes in one do not influence the other.

Multiple Regression

After getting a general sense of the data, we conduct multiple regression analysis to examine the relationship between the dependent variable, Median House Value (MEDHVAL), and the predictors, education attainment (PCTBACHMOR), number of Households Living in Poverty (NBELPOV100), percentage of Vacant Houses (PCTVACANT), and percentage of Single House Units (PCTSINGLES). Regression analysis is a statistical method to examine the relationship between a dependent variable and one or more predictors. With this type of analysis, researchers can identify the strength and direction of the relationship between variables, make predictions, and assess the significance of predictors. The model also estimates coefficients for each predictor, which represent the expected change in the dependent variable for a one-unit change in the predictor, holding other predictors constant. For this study, the multiple regression model is formulated as follows: \[ \text{LNMEDHVAL} = \beta_0 + \beta_1 \text{PCTVACANT} + \beta_2 \text{PCTSINGLES} + \beta_3 \text{PCTBACHMOR} + \beta_4 \text{log(NBELPOV100)} + \epsilon \] where LNMEDHVAL is the log-transformed median house value, PCTVACANT is the proportion of vacant housing units , PCTSINGLES is the proportion of the single family housing , PCTBACHMOR is the percentage of the residents holding bachelor’s degree or higher, and log(NBELPOV100) is the log-transformed number of households living below the poverty line.

\(\beta_0\) is the intercept, \(\beta_1\), \(\beta_2\), \(\beta_3\), and \(\beta_4\) are the coefficients for each predictor, and \(\epsilon\) is the error term. The coefficient \(\beta_1\), \(\beta_2\), \(\beta_3\), \(\beta_4\) represent the change in the log-transformed median house value for a one-unit change in the corresponding predictor, holding other predictors constant. The error term \(\epsilon\) accounts for the variability in the dependent variable that is not explained by the predictors.

Regression Assumptions

There are several assumptions associated with regression analysis that need to be met for the results to be valid. These assumptions including linearity, independence of observations, homoscedasticity, normality of residuals, no multicollinearity, and no fewer than 10 observations per predictors.

First, Linearity assumes that the relationship between the dependent variable and the predictors is linear. To verify this assumption, we made scatter plots of the dependent variable against each predictor. If the relationship appears to be linear, the assumptions was met.

Second, Independence of Observations assumes that the observations are independent of each other. There should be no spatial or temporal or other forms of dependence in the data.

Third, Homoscedasticity assumes that the variance of the residuals \(\epsilon\) is constant regardless of the values of each level of the predictors. To check this assumption, we made a scatter plot of the standardized residuals against the predicted values. If the residuals are evenly spread around zero, the assumption was met. Any patterns may indicate the presence of heteroscedasticity.

Fourth, Normality of Residuals assumes that the residuals are normally distributed. We examined the histogram of the standardized residuals to check if they are approximately normally distributed. If the histogram is bell-shaped, the assumption was met.

Fifth, No Multicollinearity assumes that the predictors are not highly correlated with each other. We calculated the correlation matrix of the predictors to check for multicollinearity. If the correlation coefficients are is not greater than 0.8 or less than -0.8, the assumption was met.

Finally, No Fewer than 10 Observations per Predictor assumes that there are at least 10 observations for each predictor in the model. Since there are over 1,700 observations in the dataset, this assumption was met.

Parameter Estimations

After verifing the assumptions of the regression, we start the regression analysis. Several parameters need to estimate here: - \(\beta_0\), which is the intercept - \(\beta_1, \dots, \beta_k\), which are coefficients of each independent variables - \(\sigma^2\), the variance of the error terms, which represents the variability in the dependent variable that is not explained by the predictors.

The least squares method estimates the coefficients by minimizing the sum of squared errors (SSE), which is the sum of the squared differences between the observed values and the predicted values. The formula for SSE is:

\[ \text{SSE} = \sum_{i=1}^{n} (y_i - \hat{y}_i)^2 = \sum_{i=1}^{n} (y_i - \hat{\beta}_0 - \hat{\beta}_1 x_{i1} - \dots - \hat{\beta}_k x_{ik})^2 \] where \(y_i\) is the actual or observed value, \(\hat{y}_i\) is the predicted value, \(\hat{\beta}_0, \dots, \hat{\beta}_k\) are the estimated coefficients, and \(x_{i1}, \dots, x_{ik}\) are the predictor values for the \(i\)-th observation.

With the Error Sum of Squares (SSE) calculated, we can estimate the variance of the error terms \(\sigma^2\) using the formula:

\[ \sigma^2 = \frac{\text{SSE}}{n - (k+1)} \] where \(n\) is the number of observations and \(k\) is the number of predictors in the model.

Model Evaluation

We evaluated the model’s fitness using the coefficient of determination \(R^2\) and the adjusted \(R^2\). The coefficient of determination \(R^2\) measures the proportion of variance in the dependent variable that is explained by the predictors. The adjusted \(R^2\) adjusts the \(R^2\) value based on the number of predictors, which help to provide a more accurate measure of model fit for multiple regression, as the increase the number of predictors can artificially inflate the \(R^2\) value.

To obtain \(R^2\), \(SST\) or the total sum of squares, needed to be calculated first. The \(SST\) measures the total variance in the dependent variable, given by: \[ SST = \sum_{i=1}^{n} (y_i - \bar{y})^2 \] Where \(y_i\) is the observed value, and \(\bar{y}\) is the mean of the observed value. Then, the \(R^2\) can be obtained by: \[ R^2 = 1 - \frac{SSE}{SST} \] After that, \(R^2\) is adjusted as follows based on the number of observations \(n\) and the number of predictors \(k\): \[ R^2_{\text{adj}} = 1 - \frac{(1 - R^2)(n - 1)}{n - k - 1} \]

Hypothesis Testing

In this analysis. several hypothesis tests are conducted to determine the significance of the model and its predictors. The overall significance of the model is assessed using the F-ratio. It compares the variance explained by the model to the variance not explained by the model.A higher F-ratio indicates that the model explains a significant amount of variance in the dependent variable compared to the residual variance.

The null hypothesis and alternative hypothesis for the F-ratio are stated as follows:

• The null hypothesis \(H_0\) states that all coefficients are equal to zero, meaning that the predictors do not explain the variance in the dependent variable (median house value). Stated as:

  \[ H_0: \beta_1 = \beta_2 = \beta_3 = \beta_4 = 0 \]

• The alternative hypothesis \(H_a\) states that at least one coefficient is not equal to zero. In our case, this means that at least one predictor explains the variance in the dependent variable (median house value). Stated as: \[ H_a: \text{At least one } \beta_i \neq 0 \]

We also conduct a t-test to determine the significance of each individual predictor in the model:

• The null hypothesis \(H_{0i}\) states that the coefficient for the predictor i is equal to zero, meaning that the predictor does not explain the variance in the dependent variable (median house value). Stated as:

  \[ H_{0i}: \beta_i = 0 \]

• The alternative hypothesis \(H_{ai}\) states that the coefficient for the predictor i is not equal to zero, meaning that the predictor explains the variance in the dependent variable (median house value). Stated as:

  \[ H_{ai}: \beta_i \neq 0 \]

Additional Analysis

Stepwise Regression

In addition to the multiple regression analysis, we conduct stepwise regression analysis to identify the most significant predictors of median house values in Philadelphia. Stepwise regression is a method that automatically selects the best subset of predictors for the model. It involves adding or removing predictors based on their statistical significance (P-value) and the AIC (Akaike Information Criterion). The stepwise regression analysis will help us identify the most important predictors and improve the overall model fits.

However, there are several limitations to stepwise regression. First, the final stepwise regression model is not guaranteed to be optimal in any specific sense. There may be other models that are as good as, or even better than the one selected by stepwise regression, but the procedure only yeilds a single model. Second, the stepwise regression does take into account researchers’ knowledge about the predictors. Important variables may be exclude from the model if they are not included in the initially. Moreover, This method also runs the risk of Type I and Type II errors—meaning it may include unimportant variables or exclude important ones.

Cross-Validation

To evaluate the performance of the regression model, we conduct K-fold cross-validation. Cross-validation is a technique used to evluate the performance of the predictive models. In K-fold Cross-Validation, the dataset is randomly divided into K equal-sized folds. The model is trained on K-1 folds and tested on the remaining fold. This process is repeated K times, with each fold serving as the test set exactly once. The average performance across all K folds is then calculated. Cross-validation helps to assess the model’s generalization performance and reduce the risk of overfitting.

In our analysis, we used K=5, which is a common choice for cross-validation. The root mean squared error (RMSE) is used to compare the performance of different models. The RMSE measures the difference between the predicted values and the actual values. A lower value of RMSE indicates a better predictive model.

The mean square errors is calculated as follows: \[ \text{MSE} = \frac{1}{n} \sum_{i=1}^{n} (y_i - \hat{y}_i)^2 \] And the root mean square errors (RMSE) is the square root of the mean square errors (MSE): \[ \text{RMSE} = \sqrt{\text{MSE}} \]

Software and Packages

In this analysis, we utilize the R programming language, a professional statistical software widely used for data analysis and visualization.The following libraries and packages have been employed to conduct the analysis:

• tidyverse: A collection of R packages for data manipulation and visualization.
• sf: A package for spatial data processing, especially for handling shapefiles.
• ggplot2: A popular visualization package for creating high-quality map and charts, part of the tidyverse collection.
• ggcorrplot: A package designated to visualize correlation matrices using ggplot2.
• patchwork: A package for combining multiple ggplot2 plots into a single plot.
• MASS: A package provides multiple functions and datasets for statistical analysis, including linear models.
• caret: A package for creating predictive models and conducting machine learning tasks, used in cross-validation and k-fold analysis.
• kableExtra: A package for creating tables with advanced formatting options.

Results

Exploratory Results

The table below presents the summary statistics for the key variables, an overview of the central tendency and variability. The dependent variable median house value has a mean of 66,287.73, indicating that the average median house value by census tract in Philadelphia is approximately 66,287.73 dollar. The standard deviation of 60,006.08 dollar suggests that house values exhibit significant variability across census tracts.

For the predictor variables, all have standard deviations close to or greater than the mean, indicating large differences in educational attainment, economic status, housing occupancy, and housing structure composition across Philadelphia’s census block groups.

dependent_var <- "MEDHVAL"

predictors <- c("PCTBACHMOR", "NBELPOV100", "PCTVACANT", "PCTSINGLES")

summary_stats <- data %>%
  dplyr::select(all_of(c(dependent_var, predictors))) %>%
  summarise_all(list(Mean = mean, SD = sd), na.rm = TRUE) %>%
  pivot_longer(cols = everything(), names_to = "Variable", values_to = "Value") %>%
  separate(Variable, into = c("Variable", "Stat"), sep = "_") %>%
  pivot_wider(names_from = Stat, values_from = Value)



summary_stats$Variable <- recode(summary_stats$Variable,
  "MEDHVAL" = "Median House Value",
  "NBELPOV100" = "# Households Living in Poverty",
  "PCTBACHMOR" = "% of Individuals with Bachelor’s Degrees or Higher",
  "PCTVACANT" = "% of Vacant Houses",
  "PCTSINGLES" = "% of Single House Units"
)



summary_stats <- summary_stats %>%
  mutate(
    Mean = round(Mean, 2),
    SD = round(SD, 2)
  )

summary_stats <- summary_stats %>%
  arrange(Variable == "Median House Value")

predictor_rows <- which(summary_stats$Variable != "Median House Value")
dependent_rows <- which(summary_stats$Variable == "Median House Value")

# Determine the start and end rows for each group
start_pred <- min(predictor_rows)
end_pred   <- max(predictor_rows)
start_dep  <- min(dependent_rows)
end_dep    <- max(dependent_rows)

# Create the table using kable and add extra formatting
kable(summary_stats, caption = "Summary Statistics", 
      align = c("l", "l", "l"), booktabs = TRUE, escape = FALSE ) %>%
  add_header_above(c(" " = 1, "Statistics" = 2)) %>%
  kable_styling(full_width = FALSE) %>%
  group_rows("Predictors", start_pred, end_pred) %>%
  group_rows("Dependent Variable", start_dep, end_dep)%>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = TRUE)

Summary Statistics

	Statistics
Variable	Mean	SD
Predictors
% of Individuals with Bachelor’s Degrees or Higher	16.08	17.77
# Households Living in Poverty	189.77	164.32
% of Vacant Houses	11.29	9.63
% of Single House Units	9.23	13.25
Dependent Variable
Median House Value	66287.73	60006.08

Histogram Distribution

One of the assumptions of OLS regression is that variables follow a normal distribution. To assess the distribution of key variables, their histograms are presented below. As seen in the plots, all key variables exhibit left-skewness, meaning their distributions are right-skewed (positively skewed), with a concentration of smaller values and a long right tail.

longer_version<- data %>%
  pivot_longer(cols = c("MEDHVAL", "PCTBACHMOR", "NBELPOV100", "PCTVACANT", "PCTSINGLES"),
               names_to = "Variable",
               values_to = "Value")

ggplot(longer_version,aes(x = Value)) +
  geom_histogram(aes(y = ..count..), fill = "black", alpha = 0.7) +  
  facet_wrap(~Variable, scales = "free", ncol = 3, labeller = as_labeller(c(
    "MEDHVAL" = "Median House Value",
    "PCTBACHMOR" = "% with Bachelor’s Degrees or Higher",
    "NBELPOV100" = "# Households Living in Poverty",
    "PCTVACANT" = "% of Vacant Houses",
    "PCTSINGLES" = "% of Single House Units"
  ))) +  
  labs(x = "Value", y = "Count", title = "Histograms of Dependent and Predictor Variables") +
  theme_light() +   
  theme(plot.subtitle = element_text(size = 9,face = "italic"),
        plot.title = element_text(size = 12, face = "bold"), 
        axis.text.x=element_text(size=6),
        axis.text.y=element_text(size=6), 
        axis.title=element_text(size=8))

To meet the normality assumption of OLS, log transformation is applied where necessary to improve distributional symmetry. This transformation helps normalize the distributions, making them more suitable for regression analysis. The histograms of the log-transformed variables indicate a more normal-like distribution. However, some slight skewness remains: the dependent variable (Log Median House Value) and the predictor (Log % Single House Value) are slightly left-skewed, while the predictors Log % Households in Poverty and Log % Vacant Houses are slightly right-skewed.

data <- data %>%
  mutate(
    LNMEDHVAL = log(MEDHVAL),
    LNPCTBACHMOR = log(1+PCTBACHMOR),
    LNNBELPOV100 = log(1+NBELPOV100),
    LNPCTVACANT = log(1+PCTVACANT),
    LNPCTSINGLES = log(1+PCTSINGLES)
  )

longer_version2 <- data %>%
  pivot_longer(cols = c(LNMEDHVAL, LNPCTBACHMOR ,LNNBELPOV100,LNPCTVACANT, LNPCTSINGLES),
               names_to = "Variable",
               values_to = "Value")

ggplot(longer_version2,aes(x = Value)) +
  geom_histogram(aes(y = ..count..), fill = "red", alpha = 0.7) +  
  facet_wrap(~Variable, scales = "free", ncol = 3, labeller = as_labeller(c(
    "LNMEDHVAL" = "Log Median House Value",
    "LNPCTBACHMOR" = "Log % with Bachelor’s Degree",
    "LNNBELPOV100" = "Log # Households in Poverty",
    "LNPCTVACANT" = "Log % Vacant Houses",
    "LNPCTSINGLES" = "Log % Single House Units"
  ))) +  
  labs(x = "Value", y = "Count", title = "Histograms of Dependent and log transformed Predictor Variables") +
  theme_light() +   
  theme(plot.subtitle = element_text(size = 9,face = "italic"),
        plot.title = element_text(size = 12, face = "bold"), 
        axis.text.x=element_text(size=6),
        axis.text.y=element_text(size=6), 
        axis.title=element_text(size=8))

Other regression assumptions will be assessed separately in the Regression Assumption Checks section later.

Choropleth Maps

Linearity assumes that the relationship between the dependent variable and the predictors is linear. To examine this relationship, we analyze the spatial distribution of key variables. Certain predictor maps exhibit patterns similar to the dependent variable, suggesting a strong potential relationship.

Below is the choropleth map of the dependent variable, Log Transformed Median House Value. It reveals that central city areas and Germantown/Chestnut Hill have higher median house values, while North Philadelphia has lower values.

ggplot(shape) +
  geom_sf(aes(fill = LNMEDHVAL), color = "transparent") +
  scale_fill_gradientn(colors = c("#fff0f3", "#a4133c"), 
                       name = "LNMEDHVAL", 
                       na.value = "transparent") + 
  theme(legend.text = element_text(size = 9),
        legend.title = element_text(size = 10),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        plot.subtitle = element_text(size = 9, face = "italic"),
        plot.title = element_text(size = 12, face = "bold"),
        panel.background = element_blank(),
        panel.border = element_rect(colour = "grey", fill = NA, size = 0.8)) +
  labs(title = "Log Transformed Median House Value")

To understand the spatial relationship between the dependent variable and predictors, we present four key predictor maps below as well.

shpe_longer<- shape %>%
  pivot_longer(cols = c("PCTVACANT", "PCTSINGLES", "PCTBACHMOR", "LNNBELPOV"),
               names_to = "Variable",
               values_to = "Value")
custom_titles <- c(
  PCTVACANT   = "Percent of Vacant Houses",
  PCTSINGLES  = "Percent of Single House Units",
  PCTBACHMOR  = "Percent of Bachelor's Degree or Higher",
  LNNBELPOV   = "Logged Transformed Poverty Rate"
)



plot_list <- lapply(unique(shpe_longer$Variable), function(var_name) {
  data_subset <- subset(shpe_longer, Variable == var_name)
  
  ggplot(data_subset) +
    geom_sf(aes(fill = Value), color = "transparent") +
    scale_fill_gradientn(
      colors = c("#fff0f3", "#a4133c"),
      name = var_name,
      na.value = "transparent"
    ) +
    labs(title = custom_titles[[var_name]]) +
    theme(
      legend.text = element_text(size = 8),
      legend.title = element_text(size = 10),
      legend.key.size = unit(0.3, "cm"),
      axis.text.x = element_blank(),
      axis.ticks.x = element_blank(),
      axis.text.y = element_blank(),
      axis.ticks.y = element_blank(),
      plot.subtitle = element_text(size = 9, face = "italic"),
      plot.title = element_text(size = 15, face = "bold"),
      panel.background = element_blank(),
      panel.border = element_rect(colour = "grey", fill = NA, size = 0.8)
    )
})

# Combine the plots into a grid (2 columns by 2 rows)
combined_plot <- (plot_list[[1]] + plot_list[[2]]) /
                 (plot_list[[3]] + plot_list[[4]])

combined_plot

These maps also help identify whether predictors are closely related to each other, which could cause multicollinearity problems. The Percent of Bachelor’s Degree or Higher map looks very similar to the Log Transformed Median House Value map, suggesting that areas with more college graduates tend to have higher house values. Also, the Percent of Vacant Houses and Percent of Single House Units maps share a similar pattern—North Philadelphia and University City have a high percentage of vacant houses and a low percentage of single-house units, while Germantown/Chestnut Hill and Far Northeast Philadelphia show the opposite trend. This suggests that these two predictors may be strongly related, which could lead to multicollinearity issues in the regression analysis.

Examining Multicollinearity & Correlation Matrix

To further examine multicollinearity, we create a correlation matrix between predictors. All correlation values are below 0.7, indicating that no severe multicollinearity occurs in our model.

Among them, % of Single House Units and Logged Households Living in Poverty have a correlation of -0.32, suggesting some relationship but not strong enough to cause severe multicollinearity. Looking back at the choropleth maps of these two predictors, their spatial patterns are not similar. % of Single House Units and % of Individuals with a Bachelor’s Degree or Higher have a correlation of -0.3, again indicating some relationship without severe multicollinearity. The choropleth maps of these two predictors also do not exhibit similar spatial patterns. Additionally, while the choropleth maps of Percent of Vacant Houses and Percent of Single House Units share a similar pattern, as mentioned before, their correlation value is only 0.2, suggesting a weak relationship without multicollinearity.

There is no severe multicollinearity among the predictors, supporting the assumptions of OLS regression. However, if the correlation matrix had shown high correlations (\(|r| \)), we could further confirm the severity of multicollinearity using the Variance Inflation Factor (VIF) (VIF < 5 indicates low multicollinearity and is generally acceptable, while VIF \(\geq\) 5 suggests moderate to high multicollinearity that may require attention. If VIF \(\geq\) 10, severe multicollinearity is present, necessitating corrective measures).

custom_labels <- c(
  "% of Individuals with Bachelor’s Degrees or Higher" = "PCTBACHMOR",
  "% of Vacant Houses" = "PCTVACANT",
  "% of Single House Units" = "PCTSINGLES",
  "# Households Living in Poverty" = "LNNBELPOV100"
)

predictor_vars <- data[, c("PCTVACANT", "PCTSINGLES", "PCTBACHMOR", "LNNBELPOV100")]

cor_matrix <- cor(predictor_vars, use = "complete.obs", method = "pearson")

print(cor_matrix)

##               PCTVACANT PCTSINGLES PCTBACHMOR LNNBELPOV100
## PCTVACANT     1.0000000 -0.1513734 -0.2983580    0.2495470
## PCTSINGLES   -0.1513734  1.0000000  0.1975461   -0.2905159
## PCTBACHMOR   -0.2983580  0.1975461  1.0000000   -0.3197668
## LNNBELPOV100  0.2495470 -0.2905159 -0.3197668    1.0000000

rownames(cor_matrix) <- names(custom_labels)
colnames(cor_matrix) <- names(custom_labels)


ggcorrplot(cor_matrix, 
           method = "square",   
           type = "lower",      
           lab = TRUE,       
           lab_size = 3,      
           colors = c("#d73027", "white", "#1a9850"))+
    labs(title = "Correlation Matrix for all Predictor Variables") +
    theme(plot.subtitle = element_text(size = 9, face = "italic"),
        plot.title = element_text(size = 12, face = "bold"), 
        axis.text.x = element_text(size = 7),
        axis.text.y = element_text(size = 7), 
        axis.title = element_text(size = 8))

Regression Results

Parameter Estimations

Four independent variables were included in the multiple regression model to predict the log-transformed median house value (LNMEDHVAL) in Philadelphia at census block group level, including the proportion of vacant housing units (PCTVACANT), the proportion of single-family housing units (PCTSINGLES), the proportion of residents with a bachelor’s degree or higher (PCTBACHMOR), and the number of households living below the poverty line (LNNBELPOV100). The regression results are presented below:

fit <- lm(LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + LNNBELPOV100, data=data)
summary(fit)

## 
## Call:
## lm(formula = LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + 
##     LNNBELPOV100, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.25825 -0.20391  0.03822  0.21744  2.24347 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  11.1137661  0.0465330 238.836  < 2e-16 ***
## PCTVACANT    -0.0191569  0.0009779 -19.590  < 2e-16 ***
## PCTSINGLES    0.0029769  0.0007032   4.234 2.42e-05 ***
## PCTBACHMOR    0.0209098  0.0005432  38.494  < 2e-16 ***
## LNNBELPOV100 -0.0789054  0.0084569  -9.330  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3665 on 1715 degrees of freedom
## Multiple R-squared:  0.6623, Adjusted R-squared:  0.6615 
## F-statistic: 840.9 on 4 and 1715 DF,  p-value: < 2.2e-16

The regression output tells us that the median housing values are highly significant correlated with the proportion of vacant housing units(PCTVANT), the proportion of single-family housing units(PCTSINGLES), the proportion of residents with a bachelor’s degree or higher (PCTBACHMOR), and the log-transformed number of households living below the poverty line (LNNBELPOV), with \(p-value<0.0001\) for all four predictors.

As the proportions of vacant housing units (PCTVANT) goes up by 1 unit (1%), the median house value goes down by approximately 1.915%, with holding all other three predictors constant. In addition, the respective \(p-value\) for \(\beta_1\) is less than 0.0001, indicating that if there is no relationship between PCTVANT and the dependent variable (i.e., if the null hypothesis that \(\beta_1 = 0\) is true), then the probability of getting a \(\beta_1\) coefficient estimate of -0.0191569 is less than 0.0001. We can safely reject the null hypothesis \(H_0: \beta_1 = 0\) for \(H_a: \beta_1 \neq 0\).

\[(e^{\beta} - 1) \times 100\% = (e^{-0.0191569} - 1) \times 100\% \approx -1.915\% \]

Similarly, as the proportion of single-family housing units (PCTSINGLES) goes up by 1 unit (1%), the median house value goes up by approximately 0.298%, with holding all other three predictors constant. The \(p-value\) for \(\beta_2\) is less than 0.0001, indicating that if there is no relationship between PCTSINGLES and the dependent variable (i.e., if the null hypothesis that \(\beta_2 = 0\) is true), then the probability of getting a \(\beta_2\) coefficient estimate of 0.0029769 is less than 0.0001. We can also safely reject the null hypothesis \(H_0: \beta_2 = 0\) for \(H_a: \beta_2 \neq 0\).

\[(e^{\beta} - 1) \times 100\% = (e^{-0.0029769} - 1) \times 100\% \approx 0.298\% \]

Furthermore, as the proportion of residents with a bachelor’s degree or higher (PCTBACHMOR) goes up by 1 unit (1%), the median house value goes up by approximately 2.09%, with holding all other three predictors constant. Similar to \(\beta_2\), the \(p-value\) for \(\beta_3\) is less than 0.0001, indicating that if there is no relationship between PCTBACHMOR and the dependent variable (i.e., if the null hypothesis that \(\beta_3 = 0\) is true, then the probability of getting a \(\beta_3\) coefficient estimate of 0.0209098 is less than 0.0001. We can also safely reject the null hypothesis \(H_0: \beta_3 = 0\) for \(H_a: \beta_3 \neq 0\).

\[(e^{\beta} - 1) \times 100\% = (e^{0.0209098} - 1) \times 100\% \approx 2.09\% \]

Lastly, as 1 percent increase in the number of households living below the poverty line, the median house value goes down by approximately 0.078% , with holding all other three predictors constant. The \(p-value\) for \(\beta_4\) is also less than 0.0001 indicating that if there is no relationship between number of households below poverty line and the dependent variable (i.e., if the null hypothesis that \(\beta_4 = 0\) is true), then the probability of getting a \(\beta_4\) coefficient estimate of -0.0789054 is less than 0.0001. As stated before, we can also safely reject the null hypothesis \(H_0: \beta_4 = 0\) for \(H_a: \beta_4 \neq 0\).

\[(1.01^{\beta} - 1) \times 100\% = (1.01^{-0.0789054} - 1) \times 100\% \approx -0.078\%\]

Model Fit

fit <- lm(LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + LNNBELPOV100, data=data)
summary(fit)

## 
## Call:
## lm(formula = LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + 
##     LNNBELPOV100, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.25825 -0.20391  0.03822  0.21744  2.24347 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  11.1137661  0.0465330 238.836  < 2e-16 ***
## PCTVACANT    -0.0191569  0.0009779 -19.590  < 2e-16 ***
## PCTSINGLES    0.0029769  0.0007032   4.234 2.42e-05 ***
## PCTBACHMOR    0.0209098  0.0005432  38.494  < 2e-16 ***
## LNNBELPOV100 -0.0789054  0.0084569  -9.330  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3665 on 1715 degrees of freedom
## Multiple R-squared:  0.6623, Adjusted R-squared:  0.6615 
## F-statistic: 840.9 on 4 and 1715 DF,  p-value: < 2.2e-16

First, the residuals from the model show a reasonable distribution of errors, with a median residual of 0.03822, a minimum value of -2.26, a first quartile of -0.20, and a third quartile of 0.22, and a maximum value of 2.24. This distribution shows that the model captures the variability in the dependent variable well (logged transformed median home value), although some of the residuals are larger than 2.5 standard deviations from the mean, which are the outliers.

Moreover, the model has a multiple \(R^2\) value of 0.6623 and an adjusted \(R^2\) value of 0.6615, indicating that approximately 66% of the variance in LNMEDHVAL is explained by the predictors. The \(F-statistic\) of 840.0 with a \(p-value\) less than 0.0001, indicating that the model is statistically significant.

anova_table <- anova(fit)
anova_table

## Analysis of Variance Table
## 
## Response: LNMEDHVAL
##                Df  Sum Sq Mean Sq  F value    Pr(>F)    
## PCTVACANT       1 180.392 180.392 1343.087 < 2.2e-16 ***
## PCTSINGLES      1  24.543  24.543  182.734 < 2.2e-16 ***
## PCTBACHMOR      1 235.118 235.118 1750.551 < 2.2e-16 ***
## LNNBELPOV100    1  11.692  11.692   87.054 < 2.2e-16 ***
## Residuals    1715 230.344   0.134                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The Analysis of Variance (ANOVA) table shows additional insight into the significance of each predictor in explaining the variation in the dependent variable (LNVMEDHVAL). For proportion of vacant housing units (PCTVACANT), the sum of squares is 180.392, indicating that those amount of variation in LNVMEDHVAL that can be attributed to changes in PCTVACANT. A higher sum of squares indicates a mhigher relationship between the predictor and the dependent variable. The sum of squares for the percentage of residents with a bachelor’s degree or higher (PCTBACHMOR) explains the highest portion of the variance in LNVMEDHVAL. In contrast, the sum of squares for the log-transformed number of households living below the poverty line (LNNBELPOV100) explains the least amount of variance in LNVMEDHVAL, with only 11.69.

The total sum of squares, which represents the total variance in the dependent variable (LNVMEDHVAL), is 230.34, with a mean square error of 0.134. Overall, the ANOVA table provides a comprehensive overview of the significance of each predictor in the model and the proportion of variance explained by each predictor. It confirms that all four predictors contribute significantly to explaining the variance in the dependent variable. These findings provide strong evidence of the robustness of the model and the relevance of neighbourhood characteristics in predicting housing values.

Regression Assumptions Check

In this section, we conducted variety of analysis to check if the assumption of the linear regression were met. In the earlier section, we already checked the variable distribution and multicollinearity. Here, we will check following assumption:

• Linearity
• Normality of Residuals
• Homoscedasticity

Linearity

To further examine the linear relationship between the dependent variable and the predictors, we create four scatter plots, where the x-axis represents each predictor and the y-axis represents the dependent variable, Log Transformed Median House Value.

longer<-data %>%
  pivot_longer(cols = c("PCTBACHMOR", "LNNBELPOV100", "PCTVACANT", "PCTSINGLES"),
               names_to = "Variable",
               values_to = "Value")

ggplot(longer,aes(x = Value, y = LNMEDHVAL)) +
  geom_point(color = "black", size= 0.4) +
  geom_smooth(method = "lm", color = "red", se = FALSE) + 
  facet_wrap(~ Variable, scales = "free", labeller = as_labeller(c(
    "PCTBACHMOR" = "% with Bachelor’s Degrees or Higher",
    "LNNBELPOV100" = "Logged Households Living in Poverty",
    "PCTVACANT" = "% of Vacant Houses",
    "PCTSINGLES" = "% of Single House Units"
  )))  +
  theme_light() +   
  theme(plot.subtitle = element_text(size = 9,face = "italic"),
        plot.title = element_text(size = 12, face = "bold"), 
        axis.text.x=element_text(size=6),
        axis.text.y=element_text(size=6), 
        axis.title=element_text(size=8)) +
  labs(title = "Scatter Plots of Dependent Variable vs. Predictors", 
       x = "Predictor Value", 
       y = "Log of Median House Value")

From the plots, we observe the following trends:

• Logged Households Living in Poverty shows a negative linear relationship with the dependent variable. As the percentage of households living in poverty increases, the Log Transformed Median House Value decreases. This is reflected in its correlation coefficient \(< 0\). However, the association is not strictly linear. In the graph, the points are more scattered widely, particular at lower level of poverty, where housing values vary significantly. This indicate while poverty levels have a negative impact on housing values, the effect is not uniform across all block groups. There is a substantial variation in housing values even at similar poverty levels, which suggests that other factors may also be influencing housing values.

• % with Bachelor’s Degrees or Higher exhibits a positive linear relationship with the dependent variable. When the percentage of individuals with a bachelor’s degree or higher increases, the Log Transformed Median House Value also rises, with a correlation coefficient \(> 0\). The points on the graph form a tight band along the trend line, indicating that higher education levels within a neighborhood is highly correlated with higher housing values. The consistent upward trend suggest that educational attainment is a key driver of housing prices, and a clear linear trend implies that this variable is suitable for this linear regression analysis.

• % of Single House Units has a weaker but still positive relationship with the dependent variable. While the trend is not obvious, areas with a higher percentage of single-house units tend to have higher Log Transformed Median House Values, especially when the percentage is sufficiently large. The correlation coefficient is \(> 0\). The scatter points shows significant variability in housing values at various levels of single-house units, with some block groups having high housing values despite a low percentage of single-house units. This variability suggests that while there is a positive relationship between single-house units and housing values, other factors may also be influencing housing values at the same time.

• % of Vacant Houses has a general negative linear relationship with the dependent variable. As the percentage of vacant houses increases, the Log Transformed Median House Value decreases, with a correlation coefficient \(< 0\). However, in the graph, it reveals a more complex relationship with significant variation in housing values across different vacancy rates. There are clusters of block groups with high housing values and low vacancy rates, indicating that the relationship is not strictly linear.In the middle and lower ranges of vacancy rates, the decline in housing values becomes less pronounced. This suggests that while vacant housing is associated with lower property values, the impact diminishes or becomes less predictable as the vacancy rate changes.

In summary, the scatterplot shows that while the PCTBACHMOR and LNMEDHVAL have a clear linear relationship, the other predictors have more complex non-linear relationships with the dependent variable. This suggests that the linear regression model may not fully capture the complexity of the relationships between the predictors and median house values.

Normality of Residuals

Next, we examine the normality of residuals of the regression model. The histogram of standardized residuals provides a detail view of the of the distribution of the residuals from the OLS model, showing below:

ggplot(data, aes(x = Standardized_Residuals)) +
  geom_histogram(bins = 30, fill = "black") +
  labs(title = "Histogram of Standardized Residuals", 
       x = "Standardized Residuals", 
       y = "Frequency") +
  theme_minimal() +   
  theme(plot.subtitle = element_text(size = 9,face = "italic"),
        plot.title = element_text(size = 12, face = "bold"), 
        axis.text.x=element_text(size=6),
        axis.text.y=element_text(size=6), 
        axis.title=element_text(size=8))

ggplot(data, aes(x = Fitted, y = Standardized_Residuals)) +
  geom_point(color = "black", size= 0.4) +    
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +  
  labs(
    title = "Scatter Plot of Standardized Residuals vs Fitted Values",
    x = "Predicted Values",
    y = "Standardized Residuals"
  ) +
  theme_minimal() +   
  theme(plot.subtitle = element_text(size = 9,face = "italic"),
        plot.title = element_text(size = 12, face = "bold"), 
        axis.text.x=element_text(size=6),
        axis.text.y=element_text(size=6), 
        axis.title=element_text(size=8))

join<- data %>%
  dplyr::select(POLY_ID, Standardized_Residuals)

shape <- shape %>%
  left_join(join, by = c("POLY_ID" = "POLY_ID"))

ggplot(shape)+
  geom_sf(aes(fill = Standardized_Residuals), color = "transparent") +
  scale_fill_gradientn(colors = c("#fff0f3", "#a4133c"), 
                       name = "Std Residuals", 
                       na.value = "transparent") +  # Choose a color palette, invert direction if needed
  labs(title = "Choropleth Map of Standardized Residuals") +
  theme(legend.text = element_text(size = 9),
        legend.title = element_text(size = 10),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        plot.subtitle = element_text(size = 9, face = "italic"),
        plot.title = element_text(size = 12, face = "bold"),
        panel.background = element_blank(),
        panel.border = element_rect(colour = "grey", fill = NA, size = 0.8))

Additional Analysis

stepwise_model <-  stepAIC(fit, direction = "both")

## Start:  AIC=-3448.07
## LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + LNNBELPOV100
## 
##                Df Sum of Sq    RSS     AIC
## <none>                      230.34 -3448.1
## - PCTSINGLES    1     2.407 232.75 -3432.2
## - LNNBELPOV100  1    11.692 242.04 -3364.9
## - PCTVACANT     1    51.546 281.89 -3102.7
## - PCTBACHMOR    1   199.020 429.36 -2379.0

stepwise_model$anova

## Stepwise Model Path 
## Analysis of Deviance Table
## 
## Initial Model:
## LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + LNNBELPOV100
## 
## Final Model:
## LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + LNNBELPOV100
## 
## 
##   Step Df Deviance Resid. Df Resid. Dev       AIC
## 1                       1715   230.3435 -3448.073

lm <-  trainControl(method = "cv", number = 5)

cvlm_model <- train(LNMEDHVAL ~ PCTVACANT + PCTSINGLES + PCTBACHMOR + LNNBELPOV100, data=data, method = "lm", trControl = lm)

print(cvlm_model)

## Linear Regression 
## 
## 1720 samples
##    4 predictor
## 
## No pre-processing
## Resampling: Cross-Validated (5 fold) 
## Summary of sample sizes: 1376, 1376, 1376, 1376, 1376 
## Resampling results:
## 
##   RMSE       Rsquared   MAE      
##   0.3675194  0.6638375  0.2726991
## 
## Tuning parameter 'intercept' was held constant at a value of TRUE

cvlm_model_reduced = train(LNMEDHVAL ~ PCTVACANT + MEDHHINC, data = data, method = "lm", trControl = lm)

print(cvlm_model_reduced)

## Linear Regression 
## 
## 1720 samples
##    2 predictor
## 
## No pre-processing
## Resampling: Cross-Validated (5 fold) 
## Summary of sample sizes: 1376, 1376, 1376, 1376, 1376 
## Resampling results:
## 
##   RMSE       Rsquared   MAE      
##   0.4431696  0.5091175  0.3181327
## 
## Tuning parameter 'intercept' was held constant at a value of TRUE

Discussion

Our regression analysis of Philadelphia’s census block group data clearly shows how neighborhood characteristics affect housing values. Our model explains about 66% of the variation in the log-transformed median house value (LNMEDHVAL). This indicates that factors such as the percentage of vacant houses (PCTVACANT), the percentage of single-family homes (PCTSINGLES), the percentage of residents with at least a bachelor’s degree (PCTBACHMOR), and the log-transformed number of households living below the poverty line (LNNBELPOV100) are all important in shaping housing market outcomes.

Interpreting the Results:

Need to rewrite We cannot include log-transformed value directly into the report, see the above interpret section for the interpretation of the results.

Checking the Model Assumptions:

b) Talk about the quality of the model – that is, state if this is a good model overall (e.g., R2, F-ratio test), and what other predictors that we didn’t include in our model might be associated with our dependent variable. i. Looking at the stepwise regression results, did the final model include all 4 predictors or were some dropped? What does that tell you about the quality of the model? ii. Looking at the cross-validation results, was the RMSE better for the 4 predictor model or the 2 predictor model?

To validate our model, we performed stepwise regression using the Akaike Information Criterion (AIC) and conducted 5-fold cross-validation. The stepwise regression confirmed that the selected predictors contribute to a strong model fit, while the cross-validation results showed that the full model which includes all the four predictors performs better than a reduced model that only includes the percentage of vacant houses (PCTVACANT) and median household income. These results strengthen our confidence in the predictive power and robustness of the model.

Policy Implications: A Focus on Inclusive Growth and Stability: (one paragraph only)

The findings have important implications for urban policy, particularly in ensuring that strategies to improve housing markets do not lead to displacement or increased housing burdens for vulnerable communities. The strong negative impact of vacancy rates (PCTVACANT) on housing values suggests that reducing long-term vacancies could play a role in neighborhood stabilization. However, rather than policies that solely focus on increasing property values, interventions should prioritize equitable investment like rehabilitating vacant properties for affordable housing, expanding community land trusts, and offering support for longtime residents at risk of displacement.

Similarly, the strong relationship between educational attainment (PCTBACHMOR) and housing values highlights the importance of expanding access to high-quality education and workforce development. However, increasing education levels alone does not ensure equitable housing outcomes if housing costs rise beyond affordability for existing residents. Thus, policies should be designed to balance educational investment with affordability measures, such as rental assistance programs and first-time homebuyer support for lower-income households.

The negative association between poverty levels (LNNBELPOV100) and housing values highlights the need for housing policies that actively combat economic inequality and support affordable, stable communities. Solutions such as mixed-income housing developments, income-based tax credits, and expanded tenant protections can help ensure that economic mobility does not come at the cost of displacement. The goal should be to create stable, mixed-income communities where economic growth benefits both current residents and new arrivals.

Limitations and Future Research:

Although our model provides a strong understanding of housing market dynamics, there are several limitations that future research should address. First, while our model explains a significant portion of the variation in housing values, about one-third remains unexplained. This suggests that other factors like access to transportation, historical redlining patterns, or the presence of public amenities may play important roles and should be explored in future studies.

The potential presence of spatial autocorrelation suggests that neighboring block groups may not be entirely independent, which could bias some estimates. Incorporating spatial econometric techniques, such as spatial lag models or geographically weighted regression, could improve the accuracy of future analyses.

Possibility of Ridge or LASSO Regression

Ridge and LASSO regression are used when there are highly correlated predictors or when only a subset of features is expected to be truly important. In such cases, some predictor coefficients may become excessively large, leading to overfitting in a linear model. These methods manage overfitting by adding a penalty term to the loss function, which constrains the coefficients. However, in our model, the highest correlation between predictors is 0.32, indicating no significant multicollinearity. Additionally, no predictor exhibits excessively large coefficients. Therefore, applying Ridge or LASSO regression in this case would be unnecessary.